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Abstract 

The equation of state (EOS) of a Fermi superfluid is investigated in the 
BCS-BEC crossover at zero temperature. We discuss the EOS based on 
Monte-Carlo (MC) data and asymptotic expansions and the EOS derived 
from the extended BCS (EBCS) mean-field theory. Then we introduce a 
time-dependent density functional, based on the bulk EOS and Landau's su- 
perfluid hydrodynamics with a von Weizsacker-type correction, to study the 
free expansion of the Fermi superfluid. We calculate the aspect ratio and 
the released energy of the expanding Fermi cloud showing that MC EOS and 
EBCS EOS are both compatible with the available experimental data of 6 Li 
atoms. We find that the released energy satisfies an approximate analytical 
formula that is quite accurate in the BEC regime. For an anisotropic droplet, 
our numerical simulations show an initially faster reversal of anisotropy in the 
BCS regime, later suppressed by the BEC fluid. 
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I. INTRODUCTION 



Current experiments with a Fermi gas of 6 Li or 40 K atoms in two hyperfine spin states 
operate in the regime of deep Fermi degeneracy. The experiments are concentrated across a 
Feshbach resonance, where the s-wave scattering length ap of the interatomic Fermi- Fermi 
potential varies from large negative to large positive values. In this way it has been observed a 
crossover from a Bardeen-Cooper-Schrieffer (BCS) superfluid to a Bose-Einstein condensate 
(BEC) of molecular pairs [1-3]. 

The bulk energy per particle of a two-spin attractive Fermi gas can be expressed [1-5] 
in the BCS-BEC crossover by the following equation 

m = | m , (i) 

where kp — (S^n) 1 / 3 is the Fermi wave vector, n is the number density, and f(y) is a 
universal function of the inverse interaction parameter y = (/cpa^) -1 , with a F the Fermi- 
Fermi scattering length. The full behavior of the universal function f(y) is unknown but 
one expects that in the BCS regime (y <C —1) it has the following asymptotic behavior 

as found by Yang et al. [6,7] in 1957. In this regime the system is a Fermi gas of weakly 
bound Cooper pairs where the superfluid gap energy A is exponentially small. Instead, in 
the unitarity limit (y = 0) the energy per particle is proportional to that of a non-interacting 
Fermi gas and, from Monte-Carlo (MC) results [5], one finds 

/(0) = 0.42 ± 0.02 . (3) 

Finally, in the BEC regime (y >> 1) the system is a weakly repulsive Bose gas of molecules 
of mass rriM = 2m, density iim = n/2 and interacting with au = 0.6a F (from MC results 
[5] and 4-body theory [8]). In this BEC regime one expects the asymptotic expression 

/(l/) = l^- + 0(-^), (4) 
18ira F y y b ' 2 

as found by Lee, Yang and Huang [9], again in year 1957. 



II. MONTE-CARLO VS MEAN-FIELD 



We have recently shown [10] that the unknown universal function f(y) can be modelled 
by the analytical formula 



f(y) = ai-a 2 arctan (a 3 y ^ + f ) . 

V P2 + \y\J 



(5) 



This formula has been obtained from Monte-Carlo (MC) simulations [5] and the asymptotic 
expressions. Table 1 of Ref. [10] reports the values of the interpolating value of a±, a 2 , a 3 , 
Pi and /?2- The thermodynamical formula 

d(n£(n)) h 2 k 2 F ( y , \ 

tin) = — ^— = [m - ~ 5 f (v)J • (6) 

relates the bulk chemical potential /x to the energy per particle S. We call Monte-Carlo 
equation of state (MC EOS) the equation of state /x = /i(n, a F ) obtained from Eqs. (5) and 
(6). 

Within the mean-field theory, the chemical potential \i and the gap energy A of the uni- 
form Fermi gas are instead found by solving the following extended BCS (EBCS) equations 
[11,12] 

1 2(2m) 1 / 2 



a F nh 



A 1 / 2 [°° dyy 2 {-- 1 ) , (7) 



n 



V 2^h 3 Jo { + 



By solving these two EBCS equations one obtains the chemical potential /x as a function of 
n and aF- Note that EBCS theory does not predict the correct BEC limit: the molecules 
have scattering length a M = 2a_p instead of au = 0.6a^. We call EBCS equation of state 
(EBCS EOS) the mean-field equation of state /x = /x(n, of) obtained from Eqs. (7) and (8). 
Obviously, our MC EOS is much closer than the EBCS EOS to the MC results obtained in 
Ref. [5] with a fixed node technique. 
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For completeness, we observe that within the EBCS mean-field theory the condensate 
density n of the Fermi superfluid can be written in terms of a simple formula [11,13], given 
by 

m 3/2 



In Ref. [11] we have found that the condensate fraction is exponentially small in the BCS 
regime (y <C —1) and goes to unity in the BEC regime (y >> 1). A very recent MC calculation 
[14] has confirmed this behavior but find at the unitarity limit a condensate fraction slightly 
smaller (n /(n/2) = 0.50) than the mean-field expectation (n /(n/2) = 0.66). 

III. TIME-DEPENDENT DENSITY FUNCTIONAL FOR A FERMI 

SUPERFLUID 

We propose an action functional A which depends on the superfluid order parameter 
tp(r, t) as follows 

A = J dt d 3 r S^ih ij*d t ij + ^^*V 2 ^ - U\ij\ 2 - £(H 2 )H 2 j . (10) 

The term £ is the bulk energy per particle of the system, which is a function of the number 
density n(r,t) = \ip{r, t)\ 2 . The Laplacian term ^^*V 2 ^ accounts for corrections to the 
kinetic energy due to spatial variations. In the BCS regime, where the Fermi gas is weakly 
interacting, the Laplacian term is phenomenological and it is called von Weizsdcker correc- 
tion [16]. In the BEC regime, where the gas of molecules is Bose condensed, the Laplacian 
term is due to the symmetry-breaking of the bosonic field operator and it is referred to as 
quantum pressure. Note that in the deep BEC regime our action functional reduces to the 
Gross-Pitaevskii action functional [15]. In our calculations we set the numerical coefficient 
c of the gradient term equal to unity (c = 1), to obtain the correct quantum-pressure term 
in the BEC regime. In the BCS regime, a better phenomenological choice for the parameter 
c could be c = 1/3 as suggested by Tosi et al. [17,18], or c = 1/36 as suggested by Zaremba 
and Tso [19]. 



For the initial confining trap, we consider an axially symmmetric harmonic potential 

u(r,t) = ^ k(*)V + y 2 ) + w a (t)Vl , (ii) 



where u)j(t) = UjQ(-t), with j — 1,2,3 — p,p,z and 0(t) the step function, so that, after 
the external trap is switched off at t > 0, the Fermi cloud performs a free expansion. The 
Euler-Lagrange equation for the superfluid order parameter ^(r, t) is obtained by minimizing 
the action functional A. This leads to a time-dependent nonlinear Schrodinger equation 
(TDNLSE): 

ih dti> = (-|^V 2 + U + /i(|Vf)) V • (12) 

The nonlinear term \x is the bulk chemical potential of the system given by the MC EOS or 
the EBCS EOS. As noted previously, in the deep BEC regime this TDNLSE reduces to the 
familiar Gross-Pitaevskii equation. From the TDNLSE one deduces the Landau's hydro- 
dynamics equations of superfluids at zero temperature by setting ip{j,t) = \Jn(r, t)e l5(r '*\ 
v(r, t) = ^VS(r,t), and neglecting the term (—h 2 V 2 ^/n) / (2rriy/n) , which would vanish in 
the uniform regime. These hydrodynamics equations are 

d t n + V • (nv) = , (13) 
m d t v + V (p(n) + U(r, t) + ^m« 2 J = . (14) 

These superfluid equations differ from the hydrodynamic equations of a normal fluid in 
the superfluid velocity field being irrotational, i.e. V A v = 0, so that the vorticity term 
v A (V A v) does not appear in Eq. (14). 

By using the superfluid hydrodynamics equations, the stationary state in the trap is given 
by the Thomas-Fermi profile no(r) = (jl — U(r,0)). Here Jl, the chemical potential of 
the inhomogeneous system, is fixed by the normalization condition N — J d 3 r no(r), where 
iV is the number of Fermi atoms. By imposing that the hydrodynamics equations satisfy 
the scaling solutions for the density 

^Mwrwym)'^' (15) 



and the velocity 



v(r t) = (x^& <M z^) (16) 

v(Mj [ X b 1 (ty v b 2 (ty z h(t))> [lb) 



we obtain three differential equations for the scaling variables bj(t), with j — 1,2,3 — p, p, z. 
The dynamics is well approximated by evaluating the scaling differential equations at the 
center (r = 0) of the cloud. In this case the variables bj(t) satisfy the local scaling equations 
(LSE) 

^ t (*(*)) 



w+o&r bit) = , (17) 



k=l 
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where n(t) = n (0)/ Y\ bk(t). Clearly, the LSE depend critically on the EOS p = p(n, of). 

k=l 

The TDNLSE is solved by using a finite-difference Crank-Nicolson predictor-corrector 
method, that we developed to solve the time-dependent Gross-Pitaevskii equation [20] . Ob- 
serve that imaginary-time integration of Eq. (12) by the Crank-Nicolson method generates 
the ground-state of Bose condensates in a ring and in a double-well [21] much more ac- 
curately than the steepest descent method used in the past. The simple LSE are instead 
solved by using a standardleap-frog symplectic algorithm, succesfully applied to investigate 
the order-to-chaos transition in spatially homogeneous field theories [22]. 

In Ref. [23] we have compared our time- dependent theory with the available experimental 
data. Moreover, we have compared the full TDNLSE with the LSE by using both MC EOS 
and EBCS EOS. We have found that, using the same EOS, the TDNLSE gives results always 
very close to the LSE ones. Instead, we have found some differences between MC EOS and 
EBCS EOS. Figure 1 reports the aspect ratio and the released energy of a 6 Li cloud after 
1.4 ms expansion from the trap realized at ENS-Paris [2]. In the experiment of Ref. [2] the 
free expansion of 7- 10 4 cold 6 Li atoms has been studied for different values of y = {kpap)' 1 
around the Feshbach resonance (y = 0). Unfortunately, in this experiment the thermal 
component is not negligible and thus the comparison with the zero-temperature theory is not 
fully satisfactory. Figure 1 compares the experimental data of Ref. [2] with the LSE based on 
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both MC and EBCS equation of state. This figure shows that the aspect ratio predicted by 
the two zero-temperature theories exceeds the finite-temperature experimental results. This 
is not surprising because the thermal component tends to hide the hydrodynamic expansion 
of the superflud. On the other hand, the released energy of the atomic gas is well described 
by the two zero-temperature theories, and the mean-field theory seems more accurate, also 
probably due to the thermal component. In Fig. 1 the released energy is defined as in Ref. 
[2], i.e. on the basis of the rms widths of the cloud. By energy conservation, the actual 
released energy is instead given by 

E rel = J d 3 r£[n (r)]n (r) . (18) 

It is straightforward to obtain an analytical expression for the released energy assuming a 
power-law dependence fj, = C n 1 for the chemical potential (polytropic equation of state) 
and writing £[n (r)] ~ |/x[n (r)] = f Wo ^) 7 n o(r) 7 where 7 is the effective polytropic index, 
obtained as the logarithmic derivative of the chemical potential /i [10], namely 

1{y> »an m - !/<(») • {a> 

In this way one finds the simple approximate formula 

where e_p = h 2 kp(0) 2 /(2m) is the Fermi chemical potential at the center of the trap, with 
kp(0) = (37r 2 n (0)) . Fig. 2 shows that this simple approximate formula which neglects all 
details of the initial aspect ratio produces fair semi- quantitative agreement, in particular in 
the BEC regime, with the actual released energy obtained by solving numerically Eq. (18). 

During the free expansion of the cloud the aspect ratio in the BCS regime (y <C —1) 
is measurably different from the one of the BEC regime (y 1). In Ref. [23] we have 
predicted an interesting effect: starting with the same aspect ratio of the cloud, at small 
times {tu H < 3) the aspect ratio is larger in the BCS region; at intermediate times (tujn — 4) 
the aspect ratio is enhanced close to the unitarity limit (y = 0); eventually at large times 



(tuo H > 5) the aspect ratio becomes larger in the BEC region. Here uo H = (a^a^) 1 / 3 is the 
geometric average of the trapping frequencies. This prediction is based on the numerical 
simulation of the LSE shown in Fig. 2, where we plot the aspect ratio of the expanding cloud 
as a function of the inverse interaction parameter y = 1/(/cpOf) at successive time intervals. 
At t — the aspect ratio equals the trap anisotropy A = 0.34. Of course the detailed 
sequence of deformations depends on the experimental conditions and in particular on the 
initial anisotropy, but the qualitative trend of an initially faster reversal on the BCS side, 
later suppressed by the BEC gas, is predicted for the expansion of any initially cigar-shaped 
interacting fermionic cloud. 



IV. CONCLUSIONS 

We have shown that the free expansion of a Fermi superfluid in the BCS-BEC crossover, 
that we simulate in a hydrodynamic scheme at zero temperature, reveals interesting features. 
We have found that the Monte-Carlo equation of state and the mean-field equation of state 
give similar results for the free expansion of a two-spin Fermi gas. The two theories are 
in reasonable agreement with the experimental data, which, however, are affected by the 
presence of a thermal component. Our Monte-Carlo equation of state and time-dependent 
density functonal can be used to study many other interesting properties; for instance, 
the collective oscillations of the Fermi cloud [24,25,10], the Fermi-Bose mixtures across a 
Feshbach resonance of the Fermi-Fermi scattering length, and nonlinear effects like Bose- 
Fermi solitons and shock waves. Finally, we observe that new experimental data on collective 
oscillations [26] suggest that the Monte-Carlo equation of state is more reliable than the 
mean-field equation of state. 
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( k F a F ) 

FIG. 1: Cloud of N = 7 • 10 4 6 Li atoms after 1.4 ms expansion from the trap realized at ENS- 
Paris [2] with anisotrpy A = u z /to p = 0.34. Squares: experimental data. Solid lines: numerical 
simulation with LSE and MC EOS; dashed lines: numerical simulation with LSE and EBCS EOS. 
The released energy is normalized as E Te \/(NeF)- 
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FIG. 2: Comparison of the approximate expression of Eq. (20) (dotted line) to the actual 
released energy defined in Eq. (18) for the conditions of the ENS-Paris experiment [2] (initial 
anisotrpy A = uj z /lo p = 0.34, N = 7 ■ 10 4 ) (solid line). Both calculations assume the MC EOS. The 
actual released energy based on Eq. (18) and the EBCS EOS is also reported (dashed line). 
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FIG. 3: Four successive frames of the aspect ratio of the 6 Li Fermi cloud as a function of 
y = {Ufclf) ■ At t = the Fermi cloud is cigar-shaped with a constant aspect ratio equal to 
the initial trap anisotropy A = u z /u) p = 0.34. Solid lines: numerical simulation with LSE and MC 
EOS; dashed lines: numerical simulation with LSE and EBCS EOS. 
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